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I. INTRODUCTION 



Since the early 1990's there has been an interest in understandin g an d accessing the effects of electromigration on 
kinetic instabilities of crystal steps [H-Q an d epitaxial islands [l3j|. and on morphological stability of epitaxial 
films fl4|-fl6]. Among other applications (which were developed primarily at the microscale), electromigration has 
been also used for fabrication of nanometer-sized gaps in metallic films. Such gaps are suitable for testing of the 
conductive properties of single molecules and control of their functionalities [13]- O- For instance, Ref. [l9| describes 
fabrication of nanoscale contacts by using electromigration to thin down and finally break the epitaxially grown 
ultrathin (10 ML) Ag films wetting the Si(001) substrate. The gap between contacts can be cyclically opened and 
closed by applying electromigration current at 80 K to open the gap, and enabling surface diffusion by annealing to 
the room temperature, to close it. 

Thus for this and other emerging applications at the nanoscale [13, HH it seems important to understand and 
characterize the effects of substate wetting and electromigration that are simulaneously active in the physical system. 
This paper combines these effects in a model that is based on an evolution equation(s) for the continuous profile of 
the film surface. The focus is on wetting films with isotropic surface energy, but with anisotropic adatom mobility 
d, although the model allows any combination of wetting properties and anisotropies. We factor in and discuss 
the effects on film stability and morphological evolution of the electric field that is cither parallel, or perpendicular 
to the initial planar surface of the film, and do not limit considerations to small deviations from planarity, i.e. the 
arbitrary surface slopes and even surface overhangs are permitted by the model. Models of wetting appropriate for 
continuum-level modeling of the surface diffusion-based dynamics of solid films have been developed and discussed 
extensively primarily in the context of thin film heteroepita xy [22l-[35ll . Our analysis is based on one such model, 
called the two- layer exponential model for the surface energy [24. 128l.l29j| . [32j]-[37| (which is particularly useful when 
the surface energy is anisotropic), but other models of wetting discussed in Refs. [22j]- [35[ can be used instead, and 
the results are expected to be qualitatively similar. The goal of modeling in this paper is not to match the theoretical 
results to the experiment [l9| and thus help in understanding the experiment, but rather to provide the broad analysis 
of the interplay of two effects (wetting and electromigration) that, at least to our knowledge, has not been addressed 
in prior publications. 



II. PROBLEM STATEMENT 



We consider a 2D single-crystal film of unperturbed height Hq with the ID parametric surface 
T(x(u, t),z(u, t)), where x and z are the Cartesian coordinates of a point on a surface, t is the time and u is the 
parameter along the surface. The origin of the Cartesian reference frame is on the substrate, and along the substrate 
(the x-direction) the film is assumed infinite. The z-axis is in the direction normal to the substrate and to the initial 
planar film surface. Surface marker particles will be used for computations of surface dynamics 38]. Thus x and z 
j^dj Cl| ^ re ^ reSerL ^ ^ e coorc ^ ma ^ es °f a mar ker particle, which are governed by the two coupled parabolic PDEs 

x t = Vz s = V-z u , (1) 
9 

z t = -Vx a = -V-x u . (2) 
9 

Here (and below) the subscripts t, s, u, z and x denote partial differentiation with respect to these variables, V is the 
normal velocity of the surface which incorporates the physics of the problem, and g(u, t) = s u = z u j cos 9 = y 1 x\ + z\. 
Here s is the surface arclength and the surface orientation angle, i.e. the one that the unit surface normal makes 
with the reference crystalline direction (chosen along the z-axis). If the surface slopes are bounded at all times 
(surface does not overhang), then it is more convenient to describe surface dynamics by a single evolution PDE for 
the height function h(x,t) of the film. Eqs. (p}, ([2]) can be easily reduced to such "h-equation" , which we will use for 
analysis; however, Eqs. ([J), @ will be used for most computations in this paper. Similar parametric approach was 
used recently in Ref. [41| for the computation of a hill-and-valley structure coarsening in the presence of material 
deposition (growth) and strongly anisotropic surface energy. 

Assuming that temperature is sufficiently high and surface diffusion is operative, the normal velocity is given by 

v = K M (^>« + aEoq WW)} J . ( 3 ) 

where D is the adatoms diffusivity, the atomic volume, v the adatoms surface density, kT the Boltzmann's factor, /j, 
the surface chemical potential, M(9) the anisotropic adatom mobility, Eq the applied electric field, q > the effective 
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charge of adatoms, f(9) — sin#, if the electric field is vertical, or f(6) = cos9, if the electric field is horizontal, and 
a = ±1 is used to select either stabilizing, or destabilizing action of the field for the chosen combination of the vertical 
or horizontal orientation of the field and the mobility M (9) . The two values of a correspond to two possible field 
orientations once either the horizontal, or the vertical field direction has been chosen (that is, field directed up-down, 
or left-right). The first term in V describes the high-T surface diffusion, the second term describes surface diffusion 
enabled by electromigration [42| . 

The surface chemical potential /i is assumed to have the contributions from the surface energy and the surface 
wetting interaction with the substrate/film interface: 



fi = il [jk + 7 Z cos ( 



(4) 



where k = Q 8 = g~ 3 (z uu x u — x uu z u ) is the surface mean curvature and, in the general case, j(z,9) is the height- 
and orientation-dependent, i.e. anisotropic, surface energy. (Note again that z here stands for the shortest distance 
between the substrate and a chosen point (x,z) on a film surface; this distance is the height h(x,t) of the surface if 
there is no overhangs.) In this paper we focus on the effects due to anisotropic adatom mobility, thus we will use the 
simpler isotropic model for the surface energy [22jj-[35l]: 



7 = 70*) = 7/ + (is ~ If) exp (-z/t), 



(5) 



where 7s is the (constant) energy of a substrate/gas (or vacuum) interface, I the characteristic wetting length, and 
jf the constant energy of a crystal/gas interface (that is, of the film surface). Eq. (JS|) is the interpolation between 
the two energies. In the limit of a thick film, z — ¥ 00, only the latter energy is retained in this expression (because the 
inter-molecular forces between the substrate and the surface molecules are relatively short-ranged), and in the limit 
of a film of zero thickness only the former energy is retained. Despite that Eq. (J5J) is phenomenological, it matches 
surprisingly well the experiments and the ab-initio calculations (at least for lattice- mismatched systems) [HI, [36|, H3] ■ 
Finally, we assume the (dimensionless) anisotropic adatom mobility in the form Q 



M(6) 



1 + ficos 2 [N(0 + </))] 
I + f3 cos 2 [N4>] 



(6) 



where N is the number of symmetry axes and <j> is the angle between a symmetry direction and the average surface 
orientation. /3 is a parameter determining the strength of the anisotropy. Throughout the paper we present results 
cither for j3 = (isotropic case), or for /3 = 1 and N — 4, <f> = tt/16. For the latter set of parameters values the 
graphs of the functions M(9) and M'(9) are shown in Fig. [TJ 




-0.6 -0.4 -0.2 0.2 0.4 0.6 



FIG. 1. Anisotropic mobility M(9) (solid line) and its derivative (dashed line) for = 1, N = 4 and (j> — tt/16. 



Remark 1. In the limit z —¥ 00 (where wetting effect is not operative), the present problem is most closely related 
to the one analyzed by Schimschak and Krug in Ref . [9( . The essential difference is that these authors determine the 
electric field from the solution of the potential equation with the appropriate boundary conditions on the material 
boundaries, including the moving surface [H|,|44|]. Thus their solution for the electric field is nonlocal, unlike the local 
approximation used in this paper. The local nature of the electric field explains why we did not detect traveling wave 
solutions, which are the hallmark of Refs. [9l. |44|. 

Next, we choose I as the length scale and £ 2 /D as the time scale and write the dimensionless counterparts of Eqs. 
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©, ©, where we use same notations for dimensionless variables: 

x t = [B {M(e)fi s } s + A {M(0)f(9)} s } z s = Vz s , (7) 

z t = ~Vx s , (8) 

/i = [1 + (G — 1) exp (—z)] K+ (1 — G) exp (—z) cos9, where cos 8 = z s . (9) 

Notice that the dimensionless forms of f(9) and the parametric expression for n coincide with their dimensional 
forms, and that Eq. §5§ has been substituted in Eq. (dJ). For conciseness, we keep differentiation with respect 
to the arclength, s, in the dimensionless equations, but their most transparent forms for the computations result 
when the differentiation with respect to s is replaced with the differentiation with respect to the parameter u, using 
d/ds = (l/g)d/du. In Eqs. - © B = n 2 v lf / (kT£ 2 ) is the surface diffusion parameter, A = avnE q/(kT) is the 
strength of the electric field, and G = js/lf is the ratio of substrate to film surface energies. For wetting films G > 1, 
for non- wetting films < G < 1. Notice that A may take on positive or negative values through the parameter a. 

As was mentioned above, if there is no overhangs then Eqs. (JZ])-© can be reduced to the single dimensionless 
evolution equation for the surface height h: 



h, = B 



M(h x )(l + h 2 x 



-1/2 



/"a 



A 



M(h x ) (l + h 2 x 



-1/2 



where 



[i = [I + (G - 1) exp (-h)] k+(1-G) exp (-h) (l + h 2 x ) 



-1/2 



f(h x ) 



-h xx (l + hi) 



-3/2 



(10) 



(11) 



/ (h x ) = h x , if the electric field is vertical, or / (h x ) — 1, if the electric field is horizontal. Also M (h x ) is given by Eq. 
([6|), where 9 is replaced by arctan (h x ). 

Using values: D = 1.5 X 1CT 6 cm 2 /s, = 2x 1(T 23 cm 3 , 7/ = 2 x 10 3 erg/cm 2 , v = 10 15 cm" 2 , fcT = 1.12 x 1CT 13 
erg, £ = 3 x 10~ 8 cm (0.3 nm, or 1ML), gives B = 8. Using q = 5e (where e = 5 x 10~ 10 statcoulombs is the absolute 
value of the electron charge) and the minimum applied voltage difference AV = 5 x 10~ 3 Volts acting across the 
typical distance of the order of the film height, AL — 10 nm, gives \qEo\ = 4 x 10 erg/cm, which translates to 
\A\ = 71. Applied voltage in surface electromigration experiments at the nanoscale can be as high as 1 V [19], thus 
in this study we explored the range of field strengths 71 < \A\ < 7100. 

In the following sections we analyze several representative situations. 



III. VERTICAL ELECTRIC FIELD 

A. Linear stability analysis 

The small slope approximation (\h x \ -C 1) of Eq. ([TU|) reads: 

ht = BM(Q) [{(1 - G) exp(-h) - 1} h xxx + h x (G - 1) exp(-/i) (h xx + l)] x + 

M'(0)(G - 1) [exp(-h)hl] x + AM(0)h xx - ^AM(0)h 2 x h xx + 2AM'(0)h x h xx , (12) 

where the mobility has been linearized about the flat surface h x = 0, that is, M (h x ) = M(0) + M'(0)h X) M(0) > 0. 
The last two terms are the simplest nonlinearities from the expansion of the electromigration flux that involve M (0) 
and M'(0). Without loss of generality we will take M(0) = 1 here and elsewhere in this paper. (See Fig. [TJ When 
mobility is isotropic, i.e. /3 = Q, then M — M(0) = 1, as is seen from Eq. ([6]).) Notice that the anisotropy of the 
mobility, M'(0), does not have an effect on linear stability, as it enters in the coefficients of the nonlinear terms. 
Also note that the last three terms can be written in a conservative form similar to the terms in the first line of the 
equation and to the first term in the second line (and this is how they are implemented in the code). 

Introducing small perturbation £(x, t) in Eq. (1121) by replacing h with ho + £(x,t) (where ho — Hq/£ is the 
dimensionless unperturbed film height), linearizing in £ and assuming normal modes for £, gives the perturbation 
growth rate u(k), where k is the perturbation wavenumber: 

u(k) = -B [1 + (G - 1) cxp(-ho)} fc 4 - [B(G - 1) cxp(-ho) + A] k 2 . (13) 

Remark 2. When wetting interaction is absent (thick film: ho — > oo), Eq. (1131) reduces to the standard one, 
ui(k) = —Bk 4 — Ak 2 , which reflects the stabilizing action of the surface diffusion and either stabilizing (A > 0, electric 
field is in the positive z direction), or destabilizing (A < 0, electric field is in the negative z direction) action of the 
electric field. Such film is absolutely stable when A > 0, but when A < it is long wave- unstable. 
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1. Analysis of Eq. iTTffl) 



• Wetting films (G > 1). From Eq. (|13|) one notices that wetting films are absolutely linearly stable when A > 0, 
but they are long wave-unstable when A < —B(G — 1) cxp(— ho) < (see Fig. [5]). The short-wavelength cut-off 
wavenumber, the maximum growth rate, and the wavenumber at which the latter occurs are: 

_ l -A-B(G-l)eM~ho) l [A + B(G - 1) cxp(-fop)] 2 k^_ 

c V B[l + (G-l)exp(-&o)] ' w — 4 B[i + (G-l)exp(-ho)] ' m ° x " v^' [ ' 

k c and io max are plotted in Figs. |3] and 01 The film stability decreases with increasing ho, and this trend 
saturates around ho = 10 (3 nm). That is, for the stated field strength A = — 71 films of thickness ho > 10 do 
not "feel" the stabilizing presence of the substrate, and they are as stable as the films that do not interact with 
the substrate at all. Of course, increasing field strength \A\ makes the film less stable, but increasing G makes it 
more stable, since the substrate energy provides stabilizing effect. In Figs. |3Jc) and 01c) the rate of decrease of 
k c and ui m ax with G (i.e., the rate of stabilization) increases fast with decreasing ho. For A — — 71 at ho ~ 3.5 
and G ~ 300 the entire dispersion curve uj{k) is below the fc-axis, see the dashed curves in Figs. [3jc) and@|c), 
and thus the film is absolutely linearly stable for ho <^ 3.5 and G >^ 300. By analysing the condition for 
longwave instability, A < —B(G— 1) exp(— ho) < 0, one can get an understanding of why this happens. Indeed, 




FIG. 2. Sketch of perturbation's linear growth rate u)(k) corresponding to longwave instability of the film surface. 




FIG. 3. Wetting films and vertical electric field: Instability cut-off wavenumber k c . (a) vs. ho, A — —71, G = 2; (b) vs. A, 
ho = 10 (solid line), ho — 1 (dashed line). G = 2; (c) vs. G, ho — 10 (solid line), ho = 6 (dash-dotted line), ho = 3.5 (dashed 
line). A = -71. 

this condition is equivalent to ho > In (S J 4~ 1 (1 — G)) , and for the typical values of A and B stated above and 
for moderately large G the right-hand side of the latter inequality is negative, thus the longwave instability 
occurs for any film height ho- This is similar to the situation when wetting effects are absent, see Remark 2 (but 
of course the spectrum of unstable wavenumbers and the magnitude of the maximum growth rate are different). 
There is a critical thickness hoc below which the film is absolutely linearly stable only if G — 1 ~ |j4|, i.e. when 
G is of the order of at least one hundred. 
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(a) , (b) W 




FIG. 4. Wetting films and vertical electric field: Maximum growth rate u>. Cases (a)-(c) are as in Fig. [3] 



• Non-wetting films (0 < G < I). In this case, for A < the long- wave instability occurs for any film thickness 
and any electric field strength - which is again similar to the situation when wetting effects are absent, see 
Remark 2. When A > the long- wave instability (with k c , u max and k max numerically similar to those shown 
in Eq. (jT4j) ) occurs when the film height is less than critical, ho < fiQ C = In (^4i? _1 (I — G) _1 ). When ho > h^ c 
the film is : 1 ' " ' ' ( \^ ' ' 1 ' " g. [5] It is clear that ho c increases with increasing either 
A, or G, or 

(b) 
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3000 4000 
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FIG. 5. (Color online). Non-wetting films and vertical electric field: (a) Contour plot of the critical thickness hoc- (b) Diagonal 
cross-section (from the lower left corner to the upper right corner) of the contour plot in (a). 



B. Nonlinear surface dynamics of wetting films 

Computations with the small-slope Eq. ([12")) (where M(0) = 1, M'(0) = 0) were performed for G = 2 and varying 
field strengths A that satisfy the condition for long-wave instability, A < —B(G — I) cxp(— ho) < 0. Computational 
domain was < x < X max , Xmax = 27r/fc max with periodic boundary conditions, and the initial condition was the 
small-amplitude cosine-shaped perturbation of the flat surface ho = const. All computations produced steady-state 
solutions that have the shape of a vertically stretched cosine curve with a fairly large amplitude. However, neither 
of these steady-states were confirmed when instead fully nonlinear parametric equations (JTJ) - were computed. [? ] 
Fig. [6] shows the computed morphology. Overhangs are clearly visible, and the bottom of the "pit" flattens out as it 
approaches the substrate due to increase of the repulsive force, promoting overhangs in the vicinity. 

Similarly, when Eq. (|12p was computed with the random small-amplitude initial condition on the domain < x < 
20A maa; (the mobility was again isotropic), the result was a perpetually coarsening hill-and- valley structure. This 
was again not confirmed in the computations of parametric equations. The late time morphology computed using 
the parametric equations is shown in Fig. [Tfa). The surface develops deep "pockets" whose walls overhang and 
eventually merge, resulting in tear drop-shaped voids trapped in the solid; this can be seen, for instance, on the 
interval 40 < x < 50. Fig. [7£b) shows the typical random initial condition. 
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Anisotropic mobility resulted in morphologies that are similar to ones shown in Figs. [6] and BJa), but skewed left 
or right, depending on the sign of <fi. As we pointed in Sec. IIII Al anisotropy matters in the nonlinear stage of the 
dynamics. 




FIG. 6. Wetting films, vertical electric field, isotropic mobility: One wavelength of the surface morphology on the periodic 
domain < x < \ ma x, computed using Eq. ()12p (Af'(O) = 0) starting from the small-amplitude cosine-shaped initial condition. 
A = -71, h = 10, G = 2. 




(b) 



FIG. 7. Wetting films, vertical electric field, isotropic mobility: (a) Surface morphology computed using Eqs. © - ([9]) 
(/? = in Eq. <j6j) ) on the domain < x < 20\ m ax, starting from the small-amplitude random initial condition shown in (b). 
A = -71, ho = 10, G = 2. 



IV. HORIZONTAL ELECTRIC FIELD 



A. Linear stability analysis 

In the case of the horizontal electric field the small slope approximation of Eq. (flOf reads: 
h t = B [{(1 - G) exp(-h) - 1} h xxx +h x (G- 1) exp(-/i) (h xx + 1)] -f 



M'(0)(G - 1) [exp(-h)h 2 x ] x - Ah x h xx + AM\Q)h xx - -AM'(0)h 2 x h x 



(15) 
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Again here we retained two simplest nonlinearities from the expansion of the electromigration flux term. The coefficient 
M (0) in the third-to-last term has been set equal to one. This nonlinear term and the last (nonlinear) term in Eq. 
(fTSf do not have an effect on linear stability. Unless mobility is isotropic, the second-to-last term of Eq. (fT5|) has an 
effect on linear stability. Assuming anisotropy, the growth rate is 

u(k) = —B [1 + (G- 1) exp(-fco)] k 4 ~ [B(G - 1) exp(-/i ) + AM'(Q)] k 2 . (16) 

Comparing Eqs. (I13|) and (|16p . it is clear that one can obtain stability properties in the horizontal field case by 
replacing A in Eq. (IT51) with AM'(0). Since M'(0) ~ —2.7 (see Fig. [1} for the chosen parameter values in Eq. ([B]), 
then for instability one must have A > (the electric field is in the positive x direction). Thus, for instance, the 
stability analysis of wetting films in the paragraph 1. of Sec. IIII A II translates directly into the case of the horizontal 
field simply by replacing A in the formulas by — 2.7 A, where A > 0. Figs. EJb) and@|b) are also valid in the case of 
the horizontal field if the A values along the horizontal axes are understood as values of the product —2.7 A (thus the 
absolute characteristic A values are roughly three times smaller than in the case of the vertical field). 

B. Nonlinear surface dynamics of wetting films 

1. Periodic steady-states from sinusoidal perturbations 

Computations with either the small-slope Eq. (fT5|) . or the parametric equations of the evolution of a one 

wavelength (A)[? ], small-amplitude cosine curve-shaped perturbation on a periodic domain resulted in steady-state 
profiles which have a shape of a vertically stretched cosine curve with a fairly large amplitude. The steady-state 
profiles often displayed a more sharply peaked bottom and less curved walls than the cosine curve, and the amplitude 
is significantly smaller in the (fully nonlinear) parametric case. In fact, it can be noticed from Fig. [5] that the film is 
nowhere close to dewetting the substrate for all tested field strengths, notwithstanding that the stabilizing influence 
of the substrate is minimal for the chosen parameters values (see the discussion in Sections IIII A 11 paragraph 1, and 

MM- 

Remark 3. When Eq. (fT5)) is used in a computation, the last term there is at large responsible for the saturation of 
the surface slope and the existence of the steady state. When this nonlinear term is omitted from the equation, the 
slope increases until the computation breaks down. 

In the computations of the parametric equations, unless A is larger than approximately 4A maa: , the surface evolves 
towards the steady-state mostly by vertical stretching of the initial shape. Perturbations of larger wavelength develop 
a large-amplitude, hill- and- valley type distortions which slowly coarsen into a steady-state, cosine curve-type shape. 
Fig. [8] shows the amplitude of the steady-state profile, h max — h m i n , and the height of the profile at it's lowest point, 
hmin, vs. A. It can be seen that, at least for moderate strengths of the electric field, growth of the amplitude is 
logarithmic in the vicinity of A c , and for A > A c growth is linear. The numerically determined slope of the linear 
section of A = 71 curve is 0.208, the one of the A — 710 curve is 0.217, and the one of the A = 7100 curve is 
0.209, which suggests that the slope is insensitive (or very weakly sensitive) to the field strength. The decay of h min 
with increasing wavelength mirrors the growth of the amplitude, thus the steady-state surface shape is symmetrically 
vertically stretched about the equilibrium surface position ho — 10. All computed steady-states are stable with 
respect to imposition of random small- and large-amplitude point perturbations, which we confirmed by computing 
the dynamics of such perturbed shapes. 

2. Coarsening of random initial roughness 

We employed fully nonlinear simulations based on parametric equations for determination of the coarsening laws at 
increasing electric field strength and at variable wetting strength (characterized by values of the parameters ho and 
G). Computations were performed on the domain < x < 2Q\ max with periodic boundary conditions; all runs were 
terminated after the surface evolved into a large-amplitude hill-and-valley structure with 3-5 hills. Unless wetting is 
strong, i.e. ho is small and G is large - the case that is discussed in more detail below - the slopes of the hills are 
constant 24° during coarsening, except for the short initial period. Figs. [5] -[11] are the log-log plots of the averaged, 
over ten realizations, maximum surface amplitude h max — h m i n and the averaged mean horizontal distance X between 
valleys (kinks) vs. time. 

Remark 4- We also attempted to compute coarsening dynamics resulting from the small-slope Eq. (|15|) . While 
the hill-and-valley structure does emerge and coarsens with time, the characteristic constant hill slope is nearly 90°. 



9 




FIG. 8. Wetting films, horizontal electric field, anisotropic mobility: Amplitude of the steady state vs. the perturbation 
wavelength A for different strengths of the electric field. Upper inset - zoom into the area of the graph around A = 0.5. Lower 
inset - minimum height of the steady state profile vs. A, also for different strengths of the electric field. Arrows labeled u X m ax" 
point out the X ma x value for each of the three electric field strengths. Parameters values are ho = 10, G — 2. 



This suggests that additional nonlinear terms must be retained in Eq. (| X5[) for predictive computations and raises the 
question what terms must be retained. We had not tried to obtain an adequate nonlinear small-slope model, leaving 
this agenda to future research. 

Of course, increasing the field strength results in faster coarsening, as the times needed for coarsening into a "final" 
structure are 10 3 , 10, and 1 for field strengths A = 71, 710 and 7100, respectively. These values also point out the 
decrease of the rate of change of a coarsening rate with the increasing field strength. Other than that, the coarsening 
laws are similar for the three tested field strengths. Fits to the data in the case A = 71 are shown in Figs. IHa,b). At 
small times coarsening is fast (exponential, see Fig. EJb)), then it changes to a slower power law with the exponent 
in the range ~ 0.1 ~ 0.2. (Since the values of the amplitude's logarithm are negative, we were unable to fit the 
exponential law to the data in Fig. IHta), thus we fitted the quadratic, which results in the i ei+C2 logl ° '-type law.) 

Remarkably, when ho is decreased from 10 to 6 a very different coarsening behaviors emerge, see Figs. ITOf a.b) and 
ITlT a.b). At G = 2, after the period of slow power-law coarsening at intermediate times, the coarsening accelerates to 
exponential coarsening or i ei+C2 logl o * type coarsening at late times (Figs. HUT a.b)). (For reference, the linear stability 
in the case ho — 6 is shown in Figs. [3jc) andSJc) by the dash-dotted line.) At G = 500 we did not output a sufficient 
number of data points in the beginning of the computation, but it is expected that initial coarsening is still faster 
than the power law. The amplitude starts to decrease towards the end of the simulation, and the kink-kink distance 
coarsens fast on the entire time interval. [? ] 

It seems certain that such unusual dynamics in the strongly wetting states emerges due to nonlinear " overdamping" 
of electromigration-induced faceting instability by the surface-substrate interaction force. Toward this end, in Figs. 
[T27 a.b) the surface shapes and surface slopes for ho = 10, G — 2 and ho — 6, G = 500 cases are compared at the 
time when seven hills formed on the surface. In the former case the hills are steep, have rather uniform height and 
their slopes are almost straight lines. In the latter case they are more irregular, "rounded", and the average height 
is smaller. It is reasonable to expect that in the opposite case of non-wetting films, the surface-substrate interaction 
will instead "sharpen up" the hills, i.e. increase their slopes and make the surface structures appear more spatially 
and temporally uniform. 

Finally, we remark that the discussed coarsening laws for strong wetting cases also are qualitatively different from 
the laws governing coarsening in the absence of wettin g, b ut in the presence of deposition, attachment-detachment, 
strong surface energy anisotropy, and interface kinetics [4ll. |45| . Indeed, only the coarsening exponents shown in Fig. 
[5] (week wetting) are within the same range (0.1-0.5) for nearly the whole computational time interval, as are the 
exponents computed in the cited papers. 
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FIG. 9. Wetting films, horizontal electric field, anisotropic mobility: (a) Log-log plot of the averaged (over ten realizations) 
maximum hill height; (b) Log-log plot of averaged kink-kink (valley-valley) distance. Dots: computed data, lines: fitting curves. 
A = 71, ho = 10, G = 2. 
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FIG. 10. Same as Fig. El for A = 71, ho = 6, G = 2. 



V. CONCLUSIONS 

In this paper the effects of electromigration and wetting on thin film morphologies are discussed, based on the 
continuum model of film surface dynamics. It has been shown that wetting effect modifies significantly the stability 
properties of the film and the coarsening of electromigration-induced surface roughness. Also it has been shown that 
the small-slope evolution equations that were employed in many studies of electromigration effects on surfaces, are 
often inadequate for the description of strongly-nonlinear phases of the dynamics. It is expected that the account of 
the surface energy anisotropy and the electric field non-locality (through the solution of the moving boundary value 
problem for the electric potential) will lead to uncovering of even more complicated behaviours. 
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FIG. 11. Same as Figs. [9]and[T0] for A = 71, ho = 6, G = 500. 



(a) (b) 




20 40 60 20 40 
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FIG. 12. (Color online.) Wetting films, horizontal electric field, anisotropic mobility: (a) Surface height. Black line: A = 
71, ho = 10, G = 2, red line: A = 71, ft = 6, G = 500; (b) Surface slope. Black crosses: A = 71, ho = 10, G = 2, red 
triangles: A — 71, fro = 6, G = 500. Note that the computational domains are of different length for the two presented cases, 
since the \ ma x values differ. 
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